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Abstract 

The dynamics of a bright matter wave soHton in a quasi ID Bose-Einstein 
condensate with periodically rapidly varying trap is considered. The gov- 
erning equation is derived based on averaging over fast modulations of the 
Gross-Pitaevskii (GP) equation. This equation has the form of CP equation 
with effective potential of more complicated structure than unperturbed trap. 
For the case of inverted (expulsive) quadratic trap corresponding to unsta- 
ble GP equation, the effective potential can be stable. For the bounded in 
space trap potential it is showed that the bifurcation exists, i.e., the single 
well potential bifurcates to the triple well effective potential. Stabilization of 
BEG cloud on-site state in the temporary modulated optical lattice is found, 
(analogous to the Kapitza stabilization of the pendulum). The predictions of 
the averaged GP equation are confirmed by the numerical simulations of GP 
equation with rapid perturbations. 
PACS numbers: 42.65.-k, 42.50. Ar,42.81.Dp 
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I. INTRODUCTION 



Recently bright matter wave solitons has been observed in a Bose-Einstein condensate 
|T|J^. In the experiment described in Ref. [|I| propagation in anisotropic BEC with cigar type 
geometry was considered. In Ref. the soliton is monitored by projecting the bound state 
of ~5000 atoms into expulsive harmonic potential. The soliton was observed propagating 
without changing of the form for distances is order of ~ 1mm. The expulsive harmonic trap 
used, in the experiment, corresponds to an unstable potential for the GP equation and a 
soliton can exists only as a metastable state in BEC 0,^. It is of interest to investigate 
the dynamics of a bright matter wave soliton in inhomogeneous time dependent systems, in 
particular possible stabilization of unstable dynamics or complicated dynamics of solitons 
in BEC. 

An interesting example of such a trap is an optical trap which is indeed rapidly var- 
ied in time. It is realized for instance for a condensate in dipole trap formed by a strong 
off-resonant laser field P] . A typical model is a trap formed by a harmonic potential with dif- 
ferent frequencies which is cutting at an energy V^. Numerical simulations of the condensate 
dynamics and splits in ID BEC with a positive scattering length under periodically shaken 
trap has been performed in Ref. The same form of wave equation leads to the prop- 
agation of spatial solitons in periodically modulated parabolic waveguides. This problem 
appears in investigation of intense light beams in a nonlinear waveguide with inhomogeneous 
distribution of transverse refractive index. A typical distribution can be approximated by 
a quadratic profile. The rapid variation of this profile along the longitudinal direction leads 
to a rapidly varying quadratic potential in the NLS equation. 

The general mathematical problem is to investigate localized states for the nonlinear 
wave equation with cubic nonlinearity and with rapidly varying (not small) potential of the 
form Vo = f{x)a(t/e),e ^ 1. Also it is of interest to investigate possibilities of stabilization 
of a soliton by rapid perturbations. For multi dimensional NLS equation under rapidly 
varying in space periodic potential such a problem is discussed in . 

The paper is organized as follows: 

In section 2 we describe the procedure of obtaining the averaged equation for the GP 
equation with an external potential, that is rapidly varying in time and inhomogeneous 
in space, with the form /(x) cos(f2t), where Q = 1/e. In Section 3, we considered some 
important potentials for applications, which have the forms given by f{x) = x^, 1/ cosh(fex), 
and cos(fcx). 
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II. AVERAGED GROSS-PITAEVSKII EQUATION 



The dynamics of BEC is described by the time dependent GP equation 

^^r = -— A^ + \/(r,t)^ + ^/o|V'l'^, (1) 
2m 

where m is the atom mass, g = ATch"^ ag/m, with, as the atomic scattering length, > 
corresponds to the BEC with repulsive interaction between atoms and < to the 
attractive interaction. The trap potential is given by = muj'^{y'^ + z"^) /2 + a{t){mujfx'^ /2 + 
Vi{x,t)), where Vi{x,t) is bounded potential or the optical lattice potential, a(t) describes 
the time dependence of the potential. We will specify the form of Vi later. Below we will 
consider the cigar type condensate with uj"^ ^ ujf. Within such restrictions we can look for 
the solution of (1) with form ipi^x, y, z, t) = R{y, z)U{x, t), where R satisfies the equation 

.9 9 

n . muj ,9 9, , , ^ 
-—AR + -—{y^ + z^)R = XpR, X = hu. 
2m 2 

Averaging over the transverse mode i?(i.e. multiplying by R*, 

2 mu . muj 2\ 

^0 = ^exp( —p ), 

Tin n 

and integrating over p) we obtain the quasi ID GP for U 

i-2 11 

im = -^f^- + (^Y^ + t))U + G\U\^U, (2) 
where G = g J \R\^dxdy/ J \R\'^dxdy = {2h\as\uj). In the dimensionless variables 



t = ujt/2, X = x/l, I = \ , u = \f2\aAU, 

V mu * 

we have the governing equation 

iut + Puxx + 2a\u\'^u = a{t)f{x)u. (3) 

where we will assume a periodic modulation of the trap 

a(t) = ao + ai sin{Qt) . (4) 

For example for V = mujfx'^/2,aof{x) = {uoi/uoYx'^ . Also a = ±1 are the attractive and 
repulsive two body interactions respectively. We introduce the parameter /? = ±1 to have 
possibility to use the results for optical beam propagation. For BEC system j3 = 1. 
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The field u{x,t) can be represented in the form of sum of slowly and rapidly varying 
parts U{x,t) and C,{x,t) 



u{x,t) = U{x,t) +^{x,t). (5) 

For obtaining the equation for an averaged field we will apply the asymptotic procedure 
suggested in 0, namely we will present the rapidly varying part of the field as expansion 
on Fourier series 

^ = A sin(fit) + B cos(fit) + C sm{2nt) + D cos(20i) + ... (6) 

where A,B,C,D are functions of (x,t) that are slowly varying in the scale 0(1) functions. 
By substituting the equations (^,(|^) into (|^) we obtain the next set of equations for the 
slowly varying field and the coefficients of the expansion for the rapidly varying component 

iUt + f3U^^ + 2\U\^U + U*{A^ + B^ + C^ + D^) + 
+2U{\A\^ + \B\^ + + + {BCA* - \A\^D + ACB* + \B\^D + 

+ABC* + ^ + = af{x)U + ^f{x)A (7) 

iAt - iVLB + (3A^^.^ + A\U\^A + 2U*{BC - AD) + 2U^A* + 

Q 

2U{-DA* + CB* + BC* - AD*) + {-\A\^A + ....) = aif{x)U - y/(x)D, (8) 

inA + iBt + I3B^^ + A\U\'^B + 2U'^B* + 2U*{AC + BD) + ... 

= y/(^)^^, (9) 

iCt - 2inD + pC^^ + 4\U\'^C + 2U*AB + ... = -j-f{x)B, (10) 

2iVLC + iDt + pD^^ + A\U\^D + U*{B'^ - A^) + ... = -y/(x)A. (11) 

The parameters ai, Vt are assumed ^1. As we can find from this system, the coefficients 
can be represented in the form of expansion 

C-^ + ^ D-^ + ^ (12) 
For the coefficients of expansion we have the expressions 
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ai = -iaif{x)Ut - ai{f{x)U):rx) - 2aif(x)\U\'^U + aoaif{x)U, 

h = iaif{x)U, d, = ^^t/, ci = .... (13) 

From (|^,(|T2|),(|T3|) we obtain the averaged equation for U 

lUt + /3f/.. + 2\U\^U = aof{x)U - i^^f\x)Ut 
-Ppix){f{x)UU - 2e'f{x)\U\'U + '^f{x)U. (14) 
Here e = ai/Q. 

This equation has the conserved quantity - the number of atoms 

dx{l + \^ ' )\U\^ = const. (15) 
So, it is useful to introduce the new field V by 

v = {i + !!/!M)i/2f/. (16) 

Substituting (|16D into (|1^) and keeping the terms in order we obtain the equation 



tVt + /3V,, + 2\VfV = (aofix) + jf3[U^)?)y + O(e'). (17) 

The averaged equation has the form of a modified NLS equation with a slowly varying 
potential 

Wix) = aofix) + ^/3[/.(x)]2. (18) 

This result shows that the soliton dynamics can have more complicated character than 
in the case of BEC with slowly varying parameters. We can expect here the stabilization of 
the unstable dynamics of a soliton by rapidly varying perturbation. It is a direct analogy 
with the stabilization of systems with a few degrees of freedom under rapidly varying in 
time parameters. 



III. DYNAMICS OF DIFFERENT MODELS OF TRAP POTENTIALS 

A. A quadratic potential 

As an example, we will consider the attractive BEC with a quadratic trap potential. In 
this case we have a perturbed trap potential [— ao + aism{Qt)]f{x) with f{x) = x^ and 
(3 = 1. An averaged effective potential is given by 
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(19) 



where e = ai/Q. The case when ao > corresponds to the inverted (repulsive) harmonic 
potential case which is well known to be unstable [Q. As seen that when 2e^ > |ao| the sign 
of the effective potential is reversed and the solitary wave dynamics becomes stable. The 
frequency of oscillations of the soliton center now is 

no = 2^(2e2-ao). (20) 

This result can be also obtained by using the moments method. It is useful to introduce 
two new variables: a position of a solitary solution of the governing equation (^ 

/oo 
dxx\u{x,t)\'^, (21) 
-oo 

and a field momentum 

/oo 
dx{u*{x, t)ux{x, t) — u{x, t)M*(x, t)). (22) 
-oo 

By using governing equation (^ we come to the following set of integro-differential equations 

-tPt + 2 r dx'^^^\u{x',t)\' = 0, (23) 

iXt -(3P = 0. 

In the case of a quadratic effective potential (W{x) = ax"^) the set of equations ( p3D is 
closed and we get the equation of motion for the soliton solution center, X(t), in the form 

Xtt + A(3aX = 0. (24) 

We point out that the set (^) is exact. No additional assumption was made in obtaining 
this set of equations, and it is valid for any solitary solution of the governing equation (^. 
Also, it should be noted two remarkable facts connected with the dynamics of a solitary 
solution of NLS equation in the presence of a quadratic potential (or dynamics of BEC in a 
quadratic trap): i).In this case the equations of motion (^) for the soliton center, are exact, 
ii) after averaging, the initially quadratic potential remains to be quadratic but its sign may 
be reversed. 
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Fig. [T|. Numerical simulations of the full GP equation for f{x) = with uj = 10, oq = 0.0493 
and values of the perturbation parameter ai providing reversing of the initially inverted harmonic 
potential. Two cases are presented : a) ai = 4.7124 and b) ai = 3.1415. 

The results of numerical simulation of the full Gross-Pitaevskii equation with the modulated 
in time expulsive (inverted) harmonic potential are shown in Fig. |l]. Two cases with ai = 
4.7124 and ai = 3.1415 were considered. For both cases oq = 0.0493 and Q = 10. The 
results of computation are the following : a) ai = 4.7124, Q = 1.309, (T = 4.8) Calculation 
by Eq. (|2|) gives n^ai = 1.257, (T^ai = 5.0). b) en = 3.1415, n = 0.7705, (T = 8.155) 
Calculation by Eq. (H) gives Qcai = 0.7695, {T^ai = 8.165). 

In both cases we observe the appearance of stable oscillations of the soliton position. We 
obtain good agreement between the full simulations and the averaged equation. As seen, the 
greater is the value of e, the discrepancy we found between results of calculation of Eq.(p^ 
and the ones of full numerical simulations of NLS equation. We can expect this since the ratio 
ai/Qis not small for these values of parameters. Let us estimate the values of parameters for 
the experimental situation. In the experiment the longitudinal frequency was cui = 4Hz, 
so, at the temporal modulations of this potential, with the frequency Q = 40Hz we should 
observe a stable bright matter wave soliton. 
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B. The case of a bounded potential 

Let us consider a bounded in the space trap potential for the atomic Bose-Einstein con- 
densate 1^ under fast temporal perturbations (— ao+ai sin(i7t))/(x) with /(x) = 1/ cosh(6x) 
and (3 = 1. 

After averaging over the period of fast perturbations, the effective potential W{x) has 
the form 

cosh(6x) 2 cosh(6x)'^ 

Provided that 

e'^*' 27 

-J- > y«o (26) 

this effective potential takes a triple well structure having a central minimum and two lateral 
minima. In Fig. ^ a typical form of the potential is presented, where ao = 0.05, ai = 9, u; = 
10 and b = 0.6. (As seen, the condition ( PB| ) makes one to take large ratio of parameters ai 
and Oq.) 




Fig. ^. Effective potential for f{x) = l/cosh(6x) after temporal averaging. The parameters of 
perturbation uj = 10, = 0.05 and ai = 9. Dashed line is for the case with ai = when temporal 
perturbation is turned off. 



The two lateral minima are positioned at x = ±3.96 and the central one at the point x = 0. 
The results of numerical simulations of the Gross-Pitaevskii equation are shown in Fig. 

i 




Fig. ^. Numerical simulations of the full GP equation for the bounded potential with /(x) = 
l/cosh(6x). The case a) corresponds to the finite motion of the soliton center in the central 
minimum with Xq = 0. The case b) corresponds to the motion in one of the lateral minima with 
Xq = 3.96. For both cases uj = 10, ao = 0.05 and ai = 9. 

To describe small oscillations of the soliton solution center and to obtain equations of 
motion in general case, we will have to make some approximations. As the ansatz, we 
suppose the solution to take the form u{x,t) = u{x — X{t),t). Here and after, for the sake 
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of simplicity, we also suppose the position of the minimum at Xq = and a symmetrical 
form for the soliton solution (dynamics of which we are investigating) 

|M(x,t)|2 = |u(-x,t)p. (27) 

Now we can consider the second term of the equation (|23|) . At first we rewrite the effective 
potential as a sum of antisymmetric and symmetric parts W{x) = W"^{x) + W^{x), where 
W^{x) = \{W{x) - W{~x)) and VV^ix) = l{W{x) + W{-x)). Then, substituting W{x) 
into the second term of Eq. (p3D expanding \u{x — X{t),t)\'^ in terms of X{t) and by holding 
the first two terms, 

Kx' + X,t)|^ = Kx,t)p + X^^^^l^, (28) 

we have 



/ dx' ^\ ' u{x' + X,t) ^ = 2 / dx' — -^u{x',t) ^+ 29 

J-oo ox' J-oo ox' 

7^x') d\u{x'. 
dx' dx 



(30) 



Introducing new parameters 



^2 ^ r ,^, dW'i^')dHx',t)\' ^ _ r^dx''-^\u{x'M 



u-ov dx' dx' ' uj'^ ' 

we get the following equation of motion for the center of soliton 

Xit + cj2(X + AX) =0. (31) 

Calculation of the parameters AX and uj"^ for our case (ao = 0.05, ai = 9,^7 = 10 and 
h = 0.6 with /(x) = 1/ cosh(6x)) gives 

1) central minimum: Xq^ = 0, u = 0.2539 (T = 24.74), AX = 0, X** = 0. Results of 
numerical simulation: uj = 0.251 (or T = 25.0), AX = 0. 

2) a lateral minimum: Xq^ = 3.96, u = 0.00599 (T = 104.73), AX = 0.321, X^* = 
Xo^ + AX = 4.18. Results of numerical simulation: oo = 0.0539 (or T = 116.8), X"* = 4.3. 
Here Xq^ is the position of a given minimum of the potential W(x); Xq''* is the stationary 
point (about which oscillations occur); AX is the shift of a stationary point; oo and T are 
the frequency and period of oscillations. 

From the present analysis we conclude that the theory of averaging of GP equation over 
rapidly varying perturbations works well in the case of bounded perturbation potential even 
when the expansion parameter e is sufficient large (in our case e = 0.9). The effective trap 
potential bifurcates from the one well form to the triple well structure. This can lead to the 
splitting of single attractive BEC into three parts. 
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C. A periodic potential 



Let f{x) = cos{kx). Then the averaged equation coincide with the unperturbed one but 
with the another potential: 

F{x) = -ao cos{kx) + - cos(2A;x)). (32) 

This potential has a more complicated structure and, as a result, the motion of the soliton 
center has new properties. Let us consider the motion of one soliton in such a potential. 
The single soliton solution is 

V{x,t) = 2irisech[2ri{x - Q] exp[i-^ + Ai^r]^ + ^^)t + iS], (33) 

where z = 2r]{x — (), and ( = A^t for the unperturbed problem. The effective potential 
acting on the soliton mass center is 

Wsoi = ^ cos(A;C) - , cos(2fcC). (34) 

2?7smh(2^) 2r7smh(^) 

The equation of motion is obtained via 

Qtt = —Q^- (35) 

So, for the soliton center (, we have the following equation of motion 

ne'^k'^ , ^, Tik^a . ,, ^, 

Ctt = ■ , .... sm(2A;C) + _ . , sm(A;C). (36) 
2?7smh(2;^) 2r7smh(^) 

As seen the motion of the soliton center occurs under the action of two harmonic effective 
perturbation. It represents a close analogy with the stabilization of a pendulum with the 
rapidly oscillating support point (the Kapitza pendulum). As it is well known [§], the 
oscillations of the pendulum at the point k( = tt are unstable. The rapid perturbation 
stabilize this fixed point. The stabilization condition follows from the effective potential for 
the soliton in the periodic field (0) and it is given by 

e^fc^ > acosh(^). (37) 
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Fig. ^ Numerical simulations of the full GP equation for a spatial periodic potential f{x) = 
cos(6x) with the initial soliton position at an unstable point Xq = ir/b. The case a) (dashed line) 
corresponds to infinite motion of the soliton center when ai = with the fast perturbation turned 
off. The case b) corresponds to ai = 9 when the fast perturbation stabilizes the motion of the 
soliton center. The parameters of perturbation are lo = 10, oq = 0.05, as in the case of a bounded 
potential. 

The results of numerical simulation of the full GP equation are presented in Fig.^ The 
parameters and u are the same as for the above considered bounded potential. Fig. ^a 
depicts the dynamics of a soliton center in a periodic potential — cos(A;x) with its initial 
position ( = n/k corresponding to the unstable point when ai = (the fast perturbation 
is turned off). As seen the motion is infinite. Fig. shows the case when ai = 9 (the 
fast perturbation is turned on). In this case the motion of the soliton becomes finite and it 
oscillates with the frequency Qi = 0.303 whereas the predicted value by Eq. ( PBD is 0.308. 
One can see that the agreement with the theory is remarkable. 



IV. CONCLUSION 



We have investigated the propagation of bright matter wave soliton in the Bose- Einstein 
condensate with trap potential rapidly varying in time. 

The cases of modulated in time quadratic, bounded and periodic trap potentials have 
been analyzed. For the repulsive (unstable) trap potential it is shown that there exists a 
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critical value of modulation parameters, when the matter wave soliton is stabilized. Analo- 
gous phenomenon of stabilization of unstable fixed points is found for periodic modulations 
(optical lattices) . For the bounded potential it is shown that the effective trap potential can 
bifurcate from the one well to the triple- well structure; and so, may give rise to the splitting 
of single attractive BEC into three parts. 
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FIGURES 

FIG. 1. Numerical simulations of the full GP equation for f{x) = with u = 10, ao = 0.0493 
and values of the perturbation parameter ai providing reversing of initially inverted harmonic 
potential. Two cases are presented : a) ai = 4.7124 ; and b) ai = 3.1415. 

FIG. 2. Effective potential for /(x) = 1/ cosh after temporal averaging. The parameters 
of perturbation u = 10, ckq = 0.05 and ai = 9. Dashed line is for the case with cci = when 
temporal perturbation is turned off. 

FIG. 3. Numerical simulations of the full GP equation for /(x) = 1/ cosh(6x). The case a) 
corresponds to the finite motion of the soliton center in the central minimum with Xq = 0. The 
case b) corresponds to the motion in one of the lateral minima with Xq = 3.96. For both cases 
Lo = 10, ao = 0.05 and ai = 9. 

FIG. 4. Numerical simulations of the full GP equation for a spatial periodic potential 
/(x) = cos(6x) with the initial soliton position at an unstable point Xq = n/k. The case a) 
(dashed line) corresponds to infinite motion of the soliton center when ai = and fast perturba- 
tion is turned off. The case b) corresponds to ai = 9 when the fast perturbation stabilizes the 
motion of the soliton center. The parameters of perturbation u = 10, ao = 0.05 as in the case of 
a bounded potential. 
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